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£Sj ■ In this paper we generalize and improve a method for calculating the period of a classical oscillator 

and other integrals of physical interest, which was recently developed by some of the authors. We 
derive analytical expressions that prove to be more accurate than those commonly found in the 
literature, and test the convergence of the series produced by the approach. 
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Tj" 1 INTRODUCTION 

\ There is great interest in the development of new methods for the treatment of nonlinear problems. For example, 
1 Amore and Saenz have recently considered the problem of calculating the period of a classical oscillator with high 
precision. The method that they propose proves to be quite effective and applicable with limited effort to a large 
1 spectrum of problems. In fact, some of the problems that considered in ref. Q are textbook examples, for which 
such method provides very accurate solutions; remarkably such solutions do not involve complicated transcendental 
Cl^' functions but are expressed in terms of elementary functions. Amore and collaborators have recently discussed mainly 
two strategies for solving nonlinear problems. One of them is the direct treatment of the integral that gives the period 
of oscillation or any other relevant property of the system, and the other is based on an improved Lindstedt-Poincare 

S technique 0, 0, 0, 0| ■ in both cases the authors resort to a sort of variational perturbation theory like that often used 
in quantum mechanics and other fields of theoretical physics to treat divergent series 0]. 

The purpose of the present paper is twofold: in first place we want to generalize the method of and to express 
it in a more systematic fashion; in second place we want to extend the previous analysis to consider large orders 
and discuss the convergence of the expansions that we obtain. We investigate the systematic calculation of integrals 
■ with applications in various fields of theoretical physics, such as, for example, the period of nonlinear oscillations, the 
deflection of the light by the sun, and the precession of the perihelion of a planet orbiting around the sun. We try 
to provide simple though sufficiently accurate analytical formulas, and test the convergence of the series, which give 
such formulas at low order. In doing so, we compare our results with those in recent literature. 



THE METHOD 



Many physical problems reduce to the calculation of integrals of the kind: 

f x+ dx 

1= / ax (i) 

where Q(x) has simple zeros at X- and x+ and is nonnegative in the interval X- < x < x+. Such integrals appear 
in many branches of classical mechanics as we will shortly show. In order to derive a simple analytical expression for 
the integral / we add and subtract a function Qq(x), which satisfies the same boundary conditions, and write 



V i+A (*) ' (2) 

where A(x) ee SMzgoM 

The method that we propose consists of expanding the integrand in powers of A(x) which leads to a series of the 
form 

oo 

I=J2ln, (3) 
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where 

-1/2 \ [ x+ A n {x) 



n 



^resent metho 

simple analytical solutions for a sufficiently large number of integrals in |QJ 



dx (4) 



and ( a b ) = a\/[b\(b — a)\] is a combinatorial number. Present method proves to be practical when we can obtain 



HARMONIC APPROXIMATION 

According to what was said above about the function Q(x) we can write 

Q(x) = R(x) (x — X-) (x+ — x) , (5) 
where R(x) > in X- < x < x + . A simple suitable reference function Qo{x) for many physical problems is 

Qo(x) = to 2 (x — X-) (x + — x) (6) 
where to is an adjustable parameter, so that 

A(x) = «^ . (7) 
In order to simplify the equations still further, we introduce the change of variable 

Xjl_ -(- X— X_l — X— n 

x = — h — cos (8) 

2 2 w 

where < 9 < it. We can thus rewrite the integral J5J as 



i= x - r -A= (9) 



and the terms in the expansion become 



lu \ n 

Clearly, the integral / is independent of u>, but the partial sums 

JV 



In = ~ ( 1/2 ) A™ dO . 



n=0 

will depend on that arbitrary parameter. It is therefore reasonable to require the partial sums to be independent of 
co, which leads to the principle of minimal sensitivity (PMS) 6] that states that the optimal value of U) should satisfy 

1S7 = ( 12 ) 

If we take into account the properties of the combinatorial numbers and that Q£j = — ^ (1 + A), we can easily prove 
that 

dS N 2N+1 

1 N ■ (i-O) 

OUJ LU 

According to this equation, the PMS condition is equivalent to In — which only takes place for odd values of N in 
the real field. For even values of N we may instead resort to the alternative PMS condition 8 2 Sn /dto 2 = but it is 
not necessary as we will see later on. 
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ANHARMONIC OSCILLATORS 

The periods of many anharmonic oscillators have been widely studied and therefore they are suitable benchmarks 
for new approaches. Here we consider a particle of unit mass moving in a one dimensional anharmonic potential V(x) 
and calculate the period according to the well known expression 

f x+ V2dx 

T= / -7^= . 14 

J x _ y/E ~ V(x) 

Duffing oscillator 

A widely studied example is the Duffing oscillator, which corresponds to the potential V(x) = \x 2 + ^ x A . In this 
case we choose Q{x) — E — V(x) = {A 2 — x 2 ) [5 + f {A 2 + x 2 )] , where A > is the amplitude of the oscillations, 
and Qo = ui 2 (A 2 — x 2 ). Notice that T = \/2I and that A takes a simple form: 



lt { A 2 + x 2 )-uj 2 + \ 



1 



(15) 



^ (l + cos 2 <?)-^ 2 ' ] 
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where x = Acos9 as follows from equation (JSJ with x+ = —x_ = A. Notice that the period depends only on p = pA 2 . 
The value of u) according to the PMS to first order 



/ 4 + 3p 

WPMS = V Z V 16 ) 



yields 



dj- P cos2d - (17) 

Eq. I|17|) shows that |Apms| < 1 for all values of 9, and p so that the series © converges for all values of p. 

A most interesting feature of the PMS for this model is that hn+i ( w = ^pms) = for all n = 0, 1, . . ., and for the 
same value of Upms given by eq. (|16f) . The calculation of the even terms is straightforward and yields the compact 
expression 



v/4 + 3p ^ V/ V n y\ 2 "/ 4 + 3p 

n— 



If we choose uj 2 = 1 + p in eq. (|T3|l then eq. @ gives us a well-known exact expression for the periodQ 

4 r /2 da 

1 exact — / -, , / / 9 y (1") 

Vi + P Jo yl — K sin a 
where n = 2 {i+p) • ^ e thus obtain the alternative series expansion: 

2tt g ^-l/2 x2 

Tl = 



If /x < 0, then the potential of the oscillator exhibits two barriers and the amplitude of the motion cannot be 
larger than Al = 1/^—p,. Consequently, the exact expression for the period is valid for p > — 1. Present series <|18[) 
converges uniformly for p > — 1 whereas the series in eq. I|20() does not converge for —1 < p < —2/3. 

The exact solution x{t) for the Duffing oscillator satisfies the virial theorem x 2 = x 2 + px 4 , where f(x) stands for 
the classical expectation value of f(x(t)). It is most interesting to note that the value of w in Xo(t) = Acos(uit + <j>) 
that makes xo(t) to satisfy the virial theorem for the Duffing oscillator is \^2uipms- 

In a recent paper Pclstcr et al. (PKS) calculated the leading term of the strong-coupling expansion for the 
frequency of the Duffing oscillator by means of a series produced by the Lindstedt-Poincare method with an adjustable 
harmonic frequency In the limit /i — > oo they expand the frequency for unit amplitude as 



(21) 
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If is the Nth-partial sum for the series l|18f) we obtain the coefficient 6q to order N as 



Una 7T7T 
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2 y^JV /-I 
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-1/2 
2j 



(22) 



Fig.JT] shows the logarithmic relative error calculated with this expression, with the values of b N ^ given by Pelster et 
al. U and the results of Amore et al 0, Q ■ We also show the linear fits for the first two sets of data. 

Our series has by far the best rate of convergence. In fact, altough all three series exhibit exponential convergence 



bo 



bo 



-ct-pN 



(23) 



the slope (3 of our linear fit is much greater: pp resent = ln(9) f=a 2.1972, compared to (3 PKS f=a 1.11 jgj. Moreover, our 
expression (|22|l . which is much simpler than the one of Pelster et al 0, enables us to calculate the slope (3 exactly. 

In Fig. [3 we compare the logarithmic relative error for the frequency calculated to second order with the three 
methods as a function of p. The results for the method of Pelster et al. are obtained by means of eq. (42) of their 
paper ||. In the limit p — > oo one recovers the asymptotic error of Fig. [Tl corresponding to N = 2. It is interesting to 
notice that both the present expansion and the result of Amore et al |2t |4j yield an error which is always smaller than 
the asymptotic one; on the other hand, the results of Pelster et al [8| do not follow this rule and yield a particularly 
large error in the region of p ~ 1. We clearly appreciate that our proposal is much more convenient than Kleinert's 
square-root trick 



Quadratic-sextic oscillator 

We also consider the potential V(x) = ^x 2 + ^x 6 because it has recently been treated by means of the combination 
of the methods of linearization and harmonic balance by Wu and Li [tj (WL). These authors carry out calculations 
of low order; their best approach is 

T WL = - 24?r ,p = pA 4 (24) 



'80 + 50p + V4096 + 5120p + 925p 2 
If we choose u> = J in cq. © we obtain an exact expression for the period: 



rpexact _ f dB .^g, 

VP + 3 Jo A /l + ^ (cos 2 6» + cos 4 6») 



In this case the greatest amplitude for p < satisfies p^ = pA\ = — 1. 

The application of our method is similar to that for the Duffing oscillator discussed above. We obtain: 

. 3 + p-&uj 2 + p cos(6>) 2 + p cos(6») 4 



6 LO 2 

In order to keep our equations as simple as possible we use the optimal value of lo of first order to all orders: 



which leads to 



IPAIS 



3(5p 



(26) 



upms = -. \*<) 



cos26» + cos46»] . (28) 



The expressions thus obtained by our method through order three are not as accurate as the one of Wu and Li 
(I24|l . At fourth order we have a formula that is not more complicated than cq. (|24f) and is certainly more accurate 
for all values of p: 

(4) A /2 7 r(6097185p 4 + 37821440 / o 3 + 89272320p 2 + 94371840p+ 3774873) 
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Notice, for example that 



lim J~pT e 

p — >OC 



8.413092631 , lim yfpT WL = 8.4081 , lim ^fpT {i) = 8.41292 



p — *oo 



p — >oo 



and that 



T exact (p = -0.9) = 10.93467798 , T WL (p = -0.9) = 10.62 , T (4) (p = -0.9) = 10.67 



Our method enables us to derive explicit equations at all orders. If we calculate the integrals in eq. I|10fl exactly, 
we obtain 



r = E 



4tt 



n=0 





p 




3(5^ + 8), 



3 ii 



where 



n ki n — ki 

E E E 23fcl - 

k 1= k 2 =0 k 3 =0 



n — k\ 



Ski. 



2n-2fc 2 -4fe 3 



In Fig. @ we plot the logarithmic relative error of the period for p — > oo as a function of the order of approximation. 
Once again we see that our expansion converges exponentially. Indeed, the error decreases as e~ a ~ /3N , where it is not 
hard to see that f3 = ln(5/3) « .5108. 



Other parity— invariant anharmonic oscillators 

The anharmonic oscillators with even potentials of the form 



(30) 



exhibit many features in common. For example, the period T depends on the parameter p = p,A 2K 2 and satisfies a 
small-p expansion 



T = E a ^ 

3=0 

and a large-/? expansion 

oo 

The optimal value of lu given by our method is of the form 

2 1 + up 



(31) 



(32) 



(33) 



The optimal value of k is the same to all orders for K = 2 as we have already discussed above. It changes slowly with 
N for the sextic oscillator and more rapidly for greater values of K. In order to make the discussion simpler we just 
consider the convergence of the series for the most difficult case given by 



c = lim ^fpT 

p — >oo 



(34) 



For K = 5 the value of k given by the PMS of first order does not give us a convergent series for this quantity. The 
reason is that the maximum value of |A| is greater than 1 under such condition. Values of k given by the PMS of higher 
order correct this problem. Numerical calculation suggests that the values of k given by the PMS at increasingly 
greater order approaches the value of k for which the maximum of A equals minus one times its minimum. Such 



6 



series. For concreteness we call this value of k because it balances the positive and negative values of A in such a 
way that |A| < 1. We find that is suitable for high order calculations but it is not necessarily the most convenient 
one for the derivation of simple and accurate analytical expressions like those obtained from the PMS of low order. 



Plots of 



„(2iV+l) 



o 



vs. 



1 show that the PMS always give the best approach 
N suggest that our partial sums with Kb converge 



vs. k obtained by partial sums of order 2N 

c (2N + l)_ 

to the exact value at each order. Plots of log — 

exponentially for all anharmonic oscillators H3U|) . 

For negative values of p there is periodic motion provided that p > — 1. Plots of log 
anharmonic oscillators with K = 3,4,5 show two almost straight lines, one for even N lying always lower than the 
one for odd N. In other words, the partial sums with even N and odd N converge exponentially to the actual value 
of T(p = -0.9). 



{p = -0.9) vs N for 



Quadratic cubic oscillator 

The treatment of noneven potentials is slightly different. We illustrate the procedure by means of the simplest such 
oscillator: 

n*) = y + My- (35) 

The factorization of Q(x) = E — V(x) leads to 

Q{x) = (x - x-)(x + — x)(b + bix) , (36) 

where 



i X+X— ■£+ ~l~ x — CZ7\ 

2(4 + + x 2 _) ' 1 ~ ~ 2(x 2 + + x+x- + x 2 _) 



and 



3 x-+ x+ 

M = -o^— ; — r ■ ( 38 ) 

Z X, + X+X- + X, 



A straightforward calculation based on the change of variable (JSJl shows that the PMS at first order yields 



f— x z , — <ix+X- — X 
uj = v , (39) 

4 



2\/ ' x\ + x+x- + x 2 



and 



x\ + Ax + x- +x 2 _ 

It is clear that all the perturbation corrections of odd order vanish and we obtain the same scries as in the case of 
the Duffing oscillator: 



UJ 



The potential-energy function exhibits a maximum V(xm) = l/(6/i 2 ) at xm = —^-/l 1 ] therefore, there will be 
periodic motion for all values of energy satisfying 



1 2 (x\ + x + x_ + x 2 
6/U 2 27 (x+ + x_ 



E<T^= V 1, — ,2 ■ (42) 
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In this case we obtain the exact solution by means of the integral 

T=^T - M (43) 

that yields the series l|41f) on expansion in a Taylor series about £ = 0. 

An interesting application of the main equations for the quadratic-cubic anharmonic oscillator is the precession of 
the perihelion of a planet. As shown in a previous paper by some of the authors the expression for the angular 
precession A</> that appears in most textbooks can be rewritten as 

= 2 / , - 2tt, (44) 

Jz- y / (z+ -z)(z- z-) [1 - 2GM (z + z-+ z+)\ 

where, z = 1/r, z± = l/r± and r_ < r + are the shortest (perihelia) and largest (aphelia) distances from the sun. 
Taking into account the expressions for the quartic-cubic oscillator derived above we have 



6GM 



y/i - 3GM (z+ +zl) = \/l- -j-, (45) 



where the semilatus rectum L is given by 1/L = (z + + zJ) /2. We also obtain 



GMJ a - L 



/. >. GM (z+ — Z-) 
A = £cos6>, £ = ^-t — = - , (46) 



where a = (l/z+ + 1/z-) /2 is the semimajor axis of the ellipse. We thus obtain the series 



= 2tt 



(47) 



that is a generalization to all orders of an expression of second order developed by Amore and Saenz . 

In Fig. 01 we plot the angular precession of the orbit of a planet obtained using the exact formula and the approxi- 
mations given by partial sums of our series (|47|l for values of a very close to the Schwartzchild radius. In Fig. we plot 
the logarithm of the relative error over the angular precession. We notice that our approximation is able to reproduce 
the exact result very accurately, even for very small a. However, our series does not reproduce the singular point of 
the exact integral exactly, because this singularity is due to the appearance of the third root of the cubic polynomial 
inside the integration interval. In other words, the factor function R(0) is negative for some values of 9 when a is 
smaller than the critical value a c = 97.9173. In spite of this limitation, the partial sums of our series approach the 
exact value of the integral as we add more terms. We can improve our results by splitting the integral into two parts, 
exactly at the new zero of Q{z), treating each integral in the way indicated. However, we do not deem it necessary 
to discuss such calculation here because it does not add anything new to what was already done. 

We can apply present method to anharmonic oscillators with more than two terms in the potential-energy function; 
as an example consider 



V{x) = aix 2 + a^x 3 + a^x 4 



If we write 



then we obtain 



Q(x) = E - V(x) = (x- x_) (x + - x) (b + b x x + b 2 x 2 ) 



E = —x + x_ [a 2 + a 3 (x + + x_) + 04 (x 2 + + x + x^ + 

&0 = 0,2 + a 3 { x + + x -) + a 4 ( x + + X + X_ + X 2 _} 

b\ = a 3 + a 4 (x + + x-) 

b 2 = a 4 
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and 

_ b + bix + b 2 x 2 - oj 2 
~ ~ 2 ' 

It follows from the PMS condition I\ = that 



'1660 + 861 (x_ + x+) + 2b 2 (3x 2 _ + 2x + x- + 3x 2 + ) 

WPMS 

1 

The expression of the period to order zero is quite simple T as \/2Iq = \/2'k /uj but the expression to second order 
T w s/2 (Iq + I2), although simple, is rather long to be shown here. 



The simple pendulum 

The method just outlined also applies to non polynomial potentials. As a representative example we consider the 
simple pendulum, and without loss of generality we choose unit mass and length. We obtain particularly simple 
expressions if we expand the potential 

V(4>) = l-cos<£ 

in a Taylor series. For example, in the first approximation V(4>) ~ 2 /2 is harmonic and the period is independent of 
the amplitude A. The resulting textbook expression for the period is valid for small values of A. In order to test the 
accuracy of the approximations we consider the exact expression for the period 



Jo 



da 



1 - sin 2 4 sin 2 f 



If we keep one more term in the expansion V(<p) ~ 2 /2 — 4 /24 we can apply all the results for the Duffing 
oscillator with \i = —1/6. The leading term of the series (|18H yields in this case 

4\/27r 

T = 



V^A 2 



Taking into account terms of higher order of the series (|18|) does not improve results because the mayor source of 
error is the expansion of the potential only through fourth order. 

The next approximation in the expansion of the potential is V{<j>) « 2 /2 — 4 /24 + </> 6 /720. Present method yields 



V6 (A 4 - 24A 2 + 192) 

oj = — i '- (48) 

48 V ' 

at first order of PMS. Keeping just the leading term of the series (0 we obtain 

T= , (49) 
VA 4 - 2AA 2 + 192 

This expression gives better results than the preceding one, and in this case the term of second order improves the 
accuracy (the term of first order is zero because of the PMS). We thus obtain 

V3tt (253A 8 - 11904A 6 + 23328CL4 4 - 221184CL4 2 + 8847360) 



15 (A 4 ~ 2AA 2 + 192) 
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FIG. 1: Logarithmic plot of the error |6g — 60I/&0 as a function of the order N. Diamonds and crosses correspond to present 
results and those of Pelster et al. respectively. The dashed and dotted lines correspond to the linear fits of each data. The 
circles correspond to the results of Q. 



CONCLUSIONS 



In this paper we propose a technique for solving a class of integrals with applications in theoretical physics. We 
provide a systematic general method for writing those integrals in a form convenient for exact calculation and for their 
expansion in scries. We show that the resulting series are suitable for all values of the relevant physical parameters, 
from unstable equilibrium through the strong-coupling, and for all the anharmonic oscillators considered here. We 
easily derive simple low-order analytical expressions that are more accurate than those in the literature, and we have 
shown that the rate of convergence of our series is greater than those proposed by other authors. In some cases the 
greater convergence rate observed is probably due to the fact that the PMS works better on the period than on the 
frequency 

The approach that we propose in this paper is closely related to what is known as variationalperturbation method, 
which was extensively applied to quantum mechanics and other fields of theoretical physics |5j (see also additional 
references cited in 0). In particular, the method of nonlinear mappings which focuses only on the series itself, is 
more general that the standard variational perturbation theory. The nonlinear mappings can also be applied to the 
problems considered here. However, we have resorted to the technique outlined above because it has facilitated the 
analysis of the convergence of the series by focusing on properties of the integrals. 
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FIG. 2: Logarithmic plot of the error \u)2 — uj exac t\/ij exac t as a function of the order [iA 2 . The three curves correspond to the 
second order results for the present method, the LPLDE method of 0, and the method of Pelster et al. respectively. The 
horizontal lines are the asymptotic values of the errors (see Fig. 0. 
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N 

FIG. 3: Logarithmic plot of the error |6q — bo|/&o as a function of the order TV for the quadratic sextic oscillator. 




a(m) 

FIG. 4: Precession of the orbit of a planet (in seconds of degree) calculated assuming M — 1.97 x 10 30 kg, G/c 2 = 7.425 x 
10 -30 m/kg and e = 0.2506 as a function of a (i.e. the average of perihelion and aphelion) for values close to a black hole. The 
solid line is the exact result; the thin line corresponds to the asymptotic formula; the other lines correspond to the different 
orders to which our approximation has been applied. 
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FIG. 5: Logarithmic plot of the error over the precession of the perihelion as a function a. The different lines correspond to 
the different orders to which our approximation has been applied. 



